Feature Selection¶
Setup¶
In [1]:
!conda install mlxtend --yes
In [36]:
import math
import pandas as pd
import numpy as np
import pandas_profiling
import matplotlib.pyplot as plt
from mlxtend.feature_selection import SequentialFeatureSelector as SFS
from mlxtend.plotting import plot_sequential_feature_selection as plot_sfs
from sklearn.pipeline import Pipeline
from sklearn.linear_model import Lasso
from sklearn.linear_model import Ridge
from sklearn.linear_model import LinearRegression
from sklearn.svm import SVR
from sklearn.ensemble import RandomForestRegressor
from sklearn.model_selection import train_test_split
from sklearn.feature_selection import VarianceThreshold
from sklearn.feature_selection import SelectPercentile
from sklearn.feature_selection import f_regression
from sklearn.feature_selection import mutual_info_regression
from sklearn.feature_selection import SelectKBest
from sklearn.feature_selection import SelectFromModel
from sklearn.feature_selection import RFE
from sklearn.metrics import mean_absolute_error
from sklearn.metrics import mean_squared_error
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import StandardScaler
from feature_engine.categorical_encoders import OneHotCategoricalEncoder
Data¶
- AllState Kaggle Competition (https://www.kaggle.com/c/allstate-claims-severity)
In [3]:
df = pd.read_csv('allstate.csv')
df.shape
Out[3]:
In [4]:
df.dtypes
Out[4]:
In [5]:
df.head()
Out[5]:
EDA¶
- Profile a small sample
In [6]:
pandas_profiling.ProfileReport(df[0:5000].copy())
Out[6]:
Pipeline¶
In [7]:
data = df[0:5000].copy()
y = data['loss']
X = data.drop(columns=['loss', 'id'])
X = X.fillna('MISSING')
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=251)
In [8]:
catColumns = X_train.select_dtypes(include='object').columns.tolist()
numColumns = X_train.select_dtypes(include='number').columns.tolist()
In [9]:
X_train_num = X_train[numColumns]
X_test_num = X_test[numColumns]
X_train_cat = X_train[catColumns]
X_test_cat = X_test[catColumns]
scaler = StandardScaler()
ohe = OneHotCategoricalEncoder(variables=catColumns)
X_train_num = pd.DataFrame(scaler.fit_transform(X_train_num), columns = X_train_num.columns)
X_train_num = X_train_num.reset_index(drop=True)
X_test_num = pd.DataFrame(scaler.transform(X_test_num), columns = X_train_num.columns)
X_test_num = X_test_num.reset_index(drop=True)
X_train_cat = ohe.fit_transform(X_train_cat)
X_train_cat = X_train_cat.reset_index(drop=True)
X_test_cat = ohe.transform(X_test_cat)
X_test_cat = X_test_cat.reset_index(drop=True)
X_train = pd.concat([X_train_num, X_train_cat], axis=1)
X_test = pd.concat([X_test_num, X_test_cat], axis=1)
Loss¶
In [10]:
def evaluate(truth, pred):
mae = mean_absolute_error(truth, pred)
mse = mean_squared_error(truth, pred)
rmse = math.sqrt(mse)
return (mae, mse, rmse)
resultsDF = pd.DataFrame(columns=['Technique', 'Method', 'MAE', 'RMSE', 'Features'])
Baseline¶
- For each technique we fit a RandomForest Regressor
- Capture MAE, RMSE and number of features
- Number of features affect training time
In [114]:
rf = RandomForestRegressor(random_state=1337)
rf.fit(X_train, y_train)
pred = rf.predict(X_test)
mae, mse, rmse = evaluate(y_test, pred)
print('MAE : {:.2f}'.format(mae))
print('RMSE : {:.2f}'.format(rmse))
resultsDF.loc['Baseline'] = ['None', 'Baseline', mae, rmse, X_train.shape[1]]
Filter¶
Highly Correlated¶
- Highly correlated features can lead to overfitting
- We can remove highly correlated variables without the loss of much information
- EDA indicated:
- cont12 and cont11 highly correlated
- cont9 and cont1 highly correlated
In [126]:
corr_matrix = X_train.corr().abs()
upper = corr_matrix.where(np.triu(np.ones(corr_matrix.shape), k=1).astype(np.bool))
features_to_remove = [column for column in upper.columns if any(upper[column] > 0.95)]
In [127]:
print(features_to_remove) # Note 0.95 > Cont9 > 0.90
X_train_new = X_train.drop(columns=features_to_remove, axis=1)
X_test_new = X_test.drop(columns=features_to_remove, axis=1)
In [136]:
# Predict & Score
rf = RandomForestRegressor(random_state=202)
rf.fit(X_train_new, y_train)
pred = rf.predict(X_test_new)
mae, mse, rmse = evaluate(y_test, pred)
print('MAE : {:.2f}'.format(mae))
print('RMSE : {:.2f}'.format(rmse))
resultsDF.loc['Filter:Correlated'] = ['Filter', 'Correlated', mae, rmse, X_train_new.shape[1]]
Select Percentile¶
- Univariate Method
- Choose your scoring method
- F-Value between label / feature
- Mutual Information
- Other selection apporaches
- Select K Best (fixed number)
- Select From Model (choose from importance or coefficients)
F Score - Regression¶
- Use a linear model for testing the individual effect of a regressor
- Scoring function:
- Correlation between regressor & target computed
- Convert to F score
- Convert to p-value
- Take top (n) or percentile of features
- For this example. pick upper 10th percentile (~80 features)
In [15]:
selector = SelectPercentile(f_regression, percentile=10)
selector.fit(X_train, y_train)
X_train_new = selector.transform(X_train)
X_test_new = selector.transform(X_test)
In [16]:
cols = X_train.columns[selector.get_support()]
pvals = selector.pvalues_[selector.get_support()]
scores = -np.log10(pvals)
scores /= scores.max()
featuresFRegDF = pd.DataFrame({'Column': cols, 'P-Value': pvals, 'Score': scores})
featuresFRegDF.sort_values('Score', ascending=False)
featuresFRegDF = featuresFRegDF.reset_index(drop=True)
featuresFRegDF['Rank'] = featuresFRegDF.index + 1
featuresFRegDF.head()
Out[16]:
In [17]:
# Predict & Score
rf = RandomForestRegressor(random_state=999)
rf.fit(X_train_new, y_train)
pred = rf.predict(X_test_new)
mae, mse, rmse = evaluate(y_test, pred)
print('MAE : {:.2f}'.format(mae))
print('RMSE : {:.2f}'.format(rmse))
resultsDF.loc['Filter:F-Score'] = ['Filter', 'F-Score', mae, rmse, X_train_new.shape[1]]
Mutual Information¶
- Mutual information measues how much a random variable tells us about another
- For example, knowing wind speed could help us guess if a flight will be delayed
- It won't tell us perfectly
- But it could make our guess better
- The degree to to which our guess is better are explained / measured with mutual information
- Note:
- Mutual information is not correlation
- Mutual information is more like information gain
- Mutual information is a distance between probability distributions
- Correlation measures the linear relationship between two random variables
- In this example, pick upper 10th percentile (about 80 features)
In [18]:
selector = SelectPercentile(mutual_info_regression, percentile=10)
selector.fit(X_train, y_train)
X_train_new = selector.transform(X_train)
X_test_new = selector.transform(X_test)
In [19]:
cols = X_train.columns[selector.get_support()]
scores = selector.scores_[selector.get_support()]
scores /= scores.max()
featuresMutDF = pd.DataFrame({'Column': cols, 'Score': scores})
featuresMutDF.sort_values('Score', ascending=False)
featuresMutDF = featuresMutDF.reset_index(drop=True)
featuresMutDF['Rank'] = featuresMutDF.index + 1
featuresMutDF.head()
Out[19]:
In [20]:
# Predict & Score
rf = RandomForestRegressor()
rf.fit(X_train_new, y_train)
pred = rf.predict(X_test_new)
mae, mse, rmse = evaluate(y_test, pred)
print('MAE : {:.2f}'.format(mae))
print('RMSE : {:.2f}'.format(rmse))
resultsDF.loc['Filter:MutualInformation'] = ['Filter', 'Mutual Information', mae, rmse, X_train_new.shape[1]]
Comparison¶
In [21]:
combineRegDF = featuresFRegDF[['Column', 'Score', 'Rank']]
combineRegDF = combineRegDF.rename(columns={'Score': 'FReg:Score', 'Rank':'FReg:Rank'})
combineMutDF = featuresMutDF[['Column', 'Score', 'Rank']]
combineMutDF = combineMutDF.rename(columns={'Score': 'Mut:Score', 'Rank':'Mut:Rank'})
featuresDF = pd.merge(left=combineRegDF, right=combineMutDF, how='outer',
left_on='Column', right_on='Column')
featuresDF
Out[21]:
Variance Threshold¶
- With variance threshold we're trying to remove features with near zero variance
- Near zero variance often equates to a poor predictor
In [120]:
selector = VarianceThreshold()
selector.fit(X_train)
selector.get_support()
X_train_new = selector.transform(X_train)
X_test_new = selector.transform(X_test)
In [121]:
# Predict & Score
rf = RandomForestRegressor(random_state=1337)
rf.fit(X_train_new, y_train)
pred = rf.predict(X_test_new)
mae, mse, rmse = evaluate(y_test, pred)
print('MAE : {:.2f}'.format(mae))
print('RMSE : {:.2f}'.format(rmse))
resultsDF.loc['Filter:VarianceThreshold'] = ['Filter', 'VarianceThreshold', mae, rmse, X_train_new.shape[1]]
Embedded¶
- Embedded method use a model to locate the best features.
- Broad approaches
- Sparsity - Lasso
- Importance - Random Forest
- Coefficients - Linear Models, SVC / SVR (different than zero)
Lasso Regression¶
- Skip fitting an intercept
In [74]:
lassoReg = Lasso(fit_intercept=False)
lassoReg.fit(X_train, y_train)
coefficients = lassoReg.coef_
In [75]:
support = coefficients > 0
cols = X_train.columns[support]
X_train_new = X_train[cols]
X_test_new = X_test[cols]
In [77]:
# Predict & Score
rf = RandomForestRegressor(random_state=999)
rf.fit(X_train_new, y_train)
pred = rf.predict(X_test_new)
mae, mse, rmse = evaluate(y_test, pred)
print('MAE : {:.2f}'.format(mae))
print('RMSE : {:.2f}'.format(rmse))
resultsDF.loc['Embedded:Lasso'] = ['Embedded', 'Lasso', mae, rmse, X_train_new.shape[1]]
Random Forest Importance¶
- Common to fit a random forest to look at importance
In [39]:
estimator = RandomForestRegressor(random_state=1337)
selector = SelectFromModel(estimator, )
selector.fit(X_train, y_train)
Out[39]:
In [40]:
cols = X_train.columns[(selector.get_support())]
X_train_new = X_train[cols]
X_test_new = X_test[cols]
In [41]:
# Predict & Score
rf = RandomForestRegressor(random_state=999)
rf.fit(X_train_new, y_train)
pred = rf.predict(X_test_new)
mae, mse, rmse = evaluate(y_test, pred)
print('MAE : {:.2f}'.format(mae))
print('RMSE : {:.2f}'.format(rmse))
resultsDF.loc['Embedded:RandomForest'] = ['Embedded', 'RandomForest', mae, rmse, X_train_new.shape[1]]
In [59]:
features = {}
for feature, importance in zip(X_train_new.columns, rf.feature_importances_):
features[feature] = importance
importanceDF = pd.DataFrame.from_dict(features, orient='index')
importanceDF.sort_values(by=0, ascending=False).head(10)
Out[59]:
In [60]:
importanceDF.sort_values(by=0, ascending=False)[0:25].plot(kind='bar', rot=90, figsize=(10, 8))
plt.title('Top 25 Important Features')
Out[60]:
Wrapper¶
Recursive Feature Elimination¶
- We start with a full model (all features)
- We then remove one feature and construct a model
- Greedily select the model with the lowest error
- Repeat until desired number of features or precisin is met
In [80]:
estimator = SVR(kernel="linear")
selector = RFE(estimator, 80, step=10, verbose=1)
selector = selector.fit(X_train, y_train)
In [81]:
X_train_new = selector.transform(X_train)
X_test_new = selector.transform(X_test)
In [82]:
# Predict & Score
rf = RandomForestRegressor(random_state=999)
rf.fit(X_train_new, y_train)
pred = rf.predict(X_test_new)
mae, mse, rmse = evaluate(y_test, pred)
print('MAE : {:.2f}'.format(mae))
print('RMSE : {:.2f}'.format(rmse))
resultsDF.loc['Wrapper:RFE'] = ['Wrapper', 'RFE', mae, rmse, X_train_new.shape[1]]
Sequential Forward Selection¶
- There isn't a forward selection in sklearn
- The mlxtend library has a robust implementation
- Parallel
- Plotting
- Target range
- Includes float (https://www.sciencedirect.com/science/article/pii/0167865594901279)
- With SFS we start with an empty model
- Build a separate model for each additional feature
- Greedily select best model
- Repeat until desired number of features or precision reached
In [106]:
estimator = LinearRegression()
sfs = SFS(estimator,
k_features=80,
forward=True,
floating=False,
verbose=2,
scoring='neg_mean_squared_error',
n_jobs=-1,
cv=0)
In [107]:
sfsModel = sfs.fit(X_train, y_train)
In [110]:
plot_sfs(sfsModel.get_metric_dict(), figsize=(10,8))
plt.title('Sequential Forward Selection')
plt.grid()
plt.show()
In [111]:
features = list(sfsModel.subsets_[80]['feature_names'])
X_train_new = X_train[features]
X_test_new = X_test[features]
In [112]:
# Predict & Score
rf = RandomForestRegressor(random_state=999)
rf.fit(X_train_new, y_train)
pred = rf.predict(X_test_new)
mae, mse, rmse = evaluate(y_test, pred)
print('MAE : {:.2f}'.format(mae))
print('RMSE : {:.2f}'.format(rmse))
resultsDF.loc['Wrapper:SFS'] = ['Wrapper', 'SFS', mae, rmse, X_train_new.shape[1]]
Sequential Forward Floating Selection¶
In [100]:
estimator = LinearRegression()
sffs = SFS(estimator,
k_features=80,
forward=True,
floating=True,
verbose=2,
scoring='neg_mean_squared_error',
n_jobs=-1,
cv=0)
In [101]:
sffsModel = sffs.fit(X_train, y_train)
In [102]:
features = list(sffsModel.subsets_[80]['feature_names'])
X_train_new = X_train[features]
X_test_new = X_test[features]
In [104]:
# Predict & Score
rf = RandomForestRegressor(random_state=999)
rf.fit(X_train_new, y_train)
pred = rf.predict(X_test_new)
mae, mse, rmse = evaluate(y_test, pred)
print('MAE : {:.2f}'.format(mae))
print('RMSE : {:.2f}'.format(rmse))
resultsDF.loc['Wrapper:SFFS'] = ['Wrapper', 'SFFS', mae, rmse, X_train_new.shape[1]]
Results¶
- RFE could have performed better (but I took steps=10)
- Note:
- These MAE values are quite close
- For a production analysis - definitely add CV or more trees to our RF evaluator
- We don't want to be unlucky on a given model
In [137]:
resultsDF.sort_values('MAE', ascending=True)
Out[137]:
Conclusion¶
- Three techniques (F.E.W.)
- Filter - Score variables and filter
- Embedded - Look at models that can identify the best features
- Wrapper - Iteratively build or remove features to find best collection
- Not PCA (But they both do a similar thing)
- Feature Selection reduces dimensionality without changing features
- PCA transforms (projects) features into a lower dimension space
- Sklearn SelectFromModel (SelectKBest, SelectPercentile) can be placed in pipelines
- No free lunch (as always)
- My default is building a full model
- Look at variable importance
- Apply RFE or Embedded (with importance)
- I'm skeptical of univariate methods in general
- Greedy
- Interactions can be meaningful
- Rather fight overfitting than underfitting
- Mlxtend worth a closer look!
In [ ]: